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ABSTRACT 

We present a study on the spatial distribution of N 2 D + and N 2 H + in thirteen 
protostellar systems. Eight of thirteen objects observed with the IRAM 30m 
telescope show relative offsets between the peak N 2 D + (J = 2 — > 1) and N 2 H + 
(J = f — > 0) emission. We highlight the case of L1157 using interferometric 
observations from the Submillimeter Array and Plateau de Bure Interferometer 
of the N 2 D+ (J = 3 — >• 2) and N 2 H+ (J = 1 ->• 0) transitions respectively. 
Depletion of N 2 D + in L1157 is clearly observed inside a radius of ~2000 AU (7") 
and the N 2 H + emission is resolved into two peaks at radii of ~1000 AU (3.5"), 
inside the depletion region of N 2 D + . Chemical models predict a depletion zone 
in N 2 D + and N 2 D + due to destruction of H 2 D + at T ~ 20 K and the evaporation 
of CO off dust grains at the same temperature. However, the abundance offsets 
of 1000 AU between the two species are not reproduced by chemical models, 
including a model that follows the infall of the protostellar envelope. The average 
abundance ratios of N 2 D + to N 2 H + have been shown to decrease as protostars 
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evolve by Emprechtinger et al., but this is the first time depletion zones of N 2 D + 
have been spatially resolved. We suggest that the difference in depletion zone 
radii for N2H + and N 2 D + is caused by either the CO evaporation temperature 
being above 20 K or an H2 ortho-to-para ratio gradient in the inner envelope. 



Introduction 



The dense (n ~10 5 cm -3 ), cold (T~10 K) gas of starless or prestellar cores present one 
of the prime conditions for the formation of deuterated molecules in the interstellar medium. 
These conditions are needed to allow the formation of H 2 D + ; reactions with this molecule and 
more highly deuterated isotopomers are the primary routes for molecules to gain molecular 
[D]/[H] ratios greater than the cosmic ratio of 10 -5 within T 



10 K gas (Watson 1974 



Langer 1985). Furthermore, at these high densities and low temperatures, CO is frozen-out 
onto dust grains, increasing the rate of deuterium fractionation in the gas-phase, as CO is 
a primary destroyer of and H 2 D + (e.g. Bergin & Tafalla 2007). A key molecular tracer 
of dense gas where CO has frozen out is N 2 H + (e.g. Casein et al. 2002). Its deuterated 



counterpart N 2 D + is found to form with abundances relative to N 2 H + of 5 - 10% in starless 



cores, with some extreme examples having 20 - 40% (Crapsi et al. 2005). 



Shortly after the gravitational collapse of a starless core begins, a Class protostar is 
formed within the dense core (Andre et al. 1993). This phase may last a few times 10 5 yr 



( Evans et al.|2009 ) and is accompanied by significant luminosity evolution and likely episodic 
accretion ( Young fc Evans||2005 Kenyon et al.|[l990 Dunham et al.pblO ). The combination 
of protostellar and accretion luminosity has a significant effect on the chemical abundances 



in the surrounding envelope (Lee et al. 2004). Regions near the protostar are heated above 
20 K (R < 1000 AU for low-mass/low-luminosity systems) and CO is evaporated off the dust 
grains back into the gas phase. Then as the protostellar luminosity increases with time, the 
region of CO evaporation increases as does the region where H 2 D + is converted back into HD 
and ELt. Thus, cold deuterium chemistry is brought to a halt in the inner envelope (Roberts 



& Millar 


2000 


Caselli et al. 


2008; 



Emprechtinger et al. (2009) showed that this process of a growing CO evaporation and 
H 2 D + destruction region will affect the N 2 D + to N 2 H + ratio and demonstrated how it could 
be used as an evolutionary indicator. The youngest Class protostars are expected to have 
higher N 2 D + to N 2 H + ratios than more evolved Class protostars. This is because the 
inside-out growth of the CO evaporation and H 2 D + destruction region causes the overall 
abundance of N 2 D + to drop relative to N 2 H + , due to the formation of N 2 D + (and other 
deuterated species) being most efficient at the highest densities. Thus, as the luminosity 
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of the protostar increases with time, the N 2 D + abundance will fall faster than the N 2 H + 
abundance. The depletion zone of N 2 D + may then be slightly larger (~200 AU) than that 
of N 2 H + , shown by the chemical models in Emprechtinger et al. ( 2009[ ). 

Depletion zones of N 2 H + have been spatially resolved previously, the most dramatic 



being IRAM 04191, with a 1600 AU radius N 2 H+ depletion zone (Belloche & Andre 2004). 



Other apparent depletion zones were pointed in Chen et al. (2007) and Tobin et al. (2011); 



depletion zones of NH3 were also seen in some sources with N 2 H + depletion regions (Tobin 



et al.||2011 ), as well as IRAM 04191 (Mangum et al. unpublished). However, most studies of 
deuterated molecules, in particular N 2 D + (e.g. Roberts k, Millar 2007), are primarily based 



on single-dish, single-point observations and do not spatially resolve or map the emission. 
N 2 D + and N 2 H + emission maps have been observed toward many starless cores and one 



protostellar core (Crapsi et al. 2004, 2005). These studies find close agreement between 



N 2 H + and N 2 D + distributions; however, a larger sample of N 2 D + maps toward protostellar 
sources are lacking. 



In conjunction with the N 2 H + kinematics survey of Tobin et al. (2011), we simultane- 
ously observed N 2 D + toward a set of thirteen low-mass protostellar systems. Furthermore, 
we execute more detailed analysis of the protostellar system L1157 where we resolve sym- 
metric depletion zones in both molecules using high-resolution interferometric observations. 
The symmetry and likely filamentary nature of the dense protostellar envelope around L1157 



e.g. Gueth et al. 2003 Looney et al. 2007; Tobin et al. 2010 Chiang et al. 2010) make it 



an ideal system for examining the detailed spatial distribution of these species. The obser- 
vations and data analysis are described in Section 2, the observational results of the spatial 
distribution of N 2 D + relative to N 2 H + in LI 157 and twelve other systems are described in 
Section 3, we compare a chemical model to L1157 in Section 4, the results are discussed and 
possible causes of the spatial offset of N 2 D + relative to N 2 H + are given in Section 5, and our 
main conclusions are given in Section 6. 



2. Observations 



This paper includes observations taken with the Submillimeter Array, the Plateau de 
Bure Interferometer and the IRAM 30m telescope. The properties of the selected sample are 



given in Table 1, with our sample being a subset of those in Tobin et al. (2011). We detail 



the observations and data reduction in the subsequent sections and list the relevant details 
of the observations in Tables 2 and 3. 
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2.1. SMA data 

L1157 was observed with the Submillimeter Array (SMA) at Mauna Kea on 2009 
September 15 in sub-compact configuration. The weather conditions were optimal with 
stable phases and less than 3 mm of precipitable water vapor, yielding a zenith opacity at 
225 GHz of 0.1; system temperatures were between 120 and 150 K throughout the track. 



The receivers were tuned to the N2D + (J = 3 — > 2) transition (z/=231. 321635 GHz; Lovas 



(1992)) in the upper-sideband and 800 MHz of each side band was allocated to continuum; 



see Table 2 for a summary of observation details. The data were edited and calibrated using 



the MIRIAD software package (Sault et al. 1995); minimal flagging of uncalibrateable data 



was necessary (i.e. phase/amplitude jumps) and the T sys correction was applied, putting 
the data on a Jansky amplitude scale. The nearby quasars 1849+670 and 1927+739 were 
used for phase and amplitude calibration, Mars was used for absolute flux calibration, and 
the quasar 3C454.3 was the bandpass calibrator. Uncertainty in absolute flux calibration is 
expected to be ~10%. 



2.2. IRAM 30m data 

Observations were taken with the IRAM 30m radio telescope on Pico Veleta using the 
EMIR receiver over several nights during 2009 October (Table 3). For L1157, we observed the 
N 2 D+ (J = 2 -»■ 1) (z/=154.217096 GHz) and (J = 3 -»■ 2) (i/=231.321966 GHz) transitions 
along with the N 2 H + (J = 3 ->■ 2) (^=279.511701 GHz) and (J = 1 -> 0) (i/=93. 173404 
GHz) transitions, see Table 3 for details. We also observed 11 other sources in the N 2 H + 
( J = 1 — > 0) and N 2 D + ( J = 2 — > 1) transitions simultaneously. System temperatures were 
between 100 K and 140 K at 93 GHz and between 110 K to 230 K at 154 GHz, see Table 3. 

All data were taken in frequency-switched On-The-Fly (OTF) mapping mode. For 
L1157 in particular, the N 2 D + (J = 3 — > 2) and N 2 H + (J = 3 — > 2) lines were mapped in 
a 90" x 90" box centered on the protostar and N 2 H+ ( J = 1 -+ 0) and N 2 D+ (J = 2 -»• 1) 
were observed in a 4' x 3' box centered 30" south of the protostar. Data for the additional 
sources were generally observed in 3' x 3' regions centered near the protostar. The OTF 
maps were constructed by fully scanning each source in both the north-south and east-west 
directions (5" steps between scan legs); the orthogonal mapping directions minimize striping 
in the final map. The maps were repeated when necessary for a higher signal-to-noise ratio. 
We conducted system temperature measurements (calibration scans) about every 10 minutes 
and about 2 hours were needed to complete each map. The pointing offset was measured 
once every two hours, preceeding the observation of each source. The offsets were generally 
~5" and this offset was stable throughout an observing session, having variations of only 
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~2" from one pointing scan to the next. The RMS pointing accuracy of the telescope is ~2" 
and our observed variations are consistent with this. The receiver offsets between the 3 mm 
and 2 mm bands were ~2" and had been measured only a month before the data were taken; 
such a small offset will not affect our results. 

The initial calibration of the OTF data to the antenna temperature scale and CLASS 
data format was performed automatically on-site by the Multichannel Imaging and calibra- 
tion software for Receiver Arrays (MIRA)[^] package. The default amplitude calibration is 
expected to be accurate within ~10%. Further data reduction was done using CLASS (part 
of GILDAS 2 ). For all molecular lines observed, the frequency switched spectra were folded 
and baseline subtracted using a second order polynomial. We then reconstructed the spectral 
map on a grid such that the FWHM of the beam was spanned by 3 pixels and each pixel is 
the average of all measurements within the FWHM of the beam, see Tables 2 and 3, as well 
Tobin et al. (2011) for additional details. Lastly, the emission maps were put on the T m b 



as 



scale where T m b = F e //T^/B e //. F e fj is the forward efficiency, accounting for losses by the 
secondary support structure, and B e jj is the main beam efficiency; values for F e ff and B e ff 
are dependent on frequency and given in Table 2. 



2.3. PdBI data 

L1157 was observed with the Plateau de Bure Interferometer in the N2H + ( J = 1 — > 0) 
transition in the D and C configurations in 2009 June and 2009 November. The data were 
edited, calibrated, and mapped using the CLIC and MAPPING components of the GILDAS 
software package. The absolute flux was derived from observations of MWC349, assuming 
a flux density of 1.15 Jy at 3 mm for that source; uncertainty in the absolute flux is ~10%. 



More details on the data reduction are given in Tobin et al. (2011) 



2.4. Single-Dish and Interferometric Combination 

As an additional step, we combined the single-dish N 2 H + data from the 30m with the 
PdBI data. We regridded the single-dish data to the phase center of the PdBI observations 
and resampled the velocity axis to correspond to the 0.125 km s _1 velocity resolution of 
the PdBI data. We did this using the MAPPING component of GILDAS, performing a 
joint deconvolution of the single-dish and interferometric data; GILDAS provides built-in 



2 http://www.iram.fr/IRAMFR/GILDAS 
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routines for combining the 30m data. The resulting image provides substantially better 



extended structure sensitivity than the default interferometric image shown in Tobin et al. 
( |20ll) . 

We also explored the possibility of combining the single-dish N2D + data with the SMA 
N 2 D + data for L1157. Unfortunately, the single-dish N 2 D + (J = 3 — > 2) data had a lower 
signal-to-noise ratio than the SMA data, prohibiting a useful combination of the data. 



2.5. Data Analysis and Hyperfine Fitting 

The nuclear spin of nitrogen causes both N2H + and N2D + to have a hyperfine emission 



line spectrum. The multiple hyperfine components and line ratios of N2H + (Caselli et al. 



1995) |Keto fc Rybicki||2010j ) and N 2 D+ ( |Dore et al7||2004[ ) enable us to simultaneously deter 



mine the line-center velocity, linewidth, optical depth, and excitation temperature at each 
point in the map. 

CLASS was used to fit the hyperfine emission line structure, which can fit multiple lines 
given the optically thin line ratios and velocity offsets. The lines were fit at each pixel in the 
spectral map using a semi-automated routine. This routine used the integrated intensity map 
(zeroth moment) of the line emission, calculated from the total line intensity in a velocity 
range and multiplied by the channel velocity width. The pixels with at least a 3a detection 
of N 2 H + or N 2 D + were selected and a CLASS script to fit the hyperfine line structure at 
each pixel was generated; the resultant fit parameters (line-center velocity, linewidth, optical 
depth, and excitation temperature) were written out to a text table. This criteria does not 
require 3a detections of all line components and there were cases in which only the strongest 
hyperfine lines are detected. Bad fits resulting from inadequate signal-to-noise were removed 
to generate the final table of pixels with good fits. The typical failure mode of the fitting is 
excessively large linewidths (>1.5 km s -1 ) and/or optical depths of 30, making the discrepant 
points easy to identify. 

The integrated intensity, excitation temperature, and optical depths from line fitting are 



then used to calculate the column densities of the data following Goldsmith & Langer (1999). 
The column density maps generally trace regions where the integrated intensity was at least 
~6cr and usually 10cr. We use the excitation temperature of the N 2 H + ( J — 1 — > 0) in the 
calculation of N 2 D + column density when the N 2 D + ( J — 2 — > 1) optical depth is less than 
1.0. This is because the excitation temperature cannot be well-constrained if the lines are 
optically thin. If the excitation temperature of N 2 D + ( J — 2 — > 1) is actually lower than that 
of N 2 H + (J — 1 — > 0), the N 2 D + abundances would systematically underestimated. When 
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the N 2 D + excitation temperatures could be constrained from the hyperfine fitting, they are 
comparable but generally lower than the values fit for N 2 H + (J — 1 — > 0). The uncertainties 
in optical depth and excitation temperature lead to column density uncertainties of order 
~ 20% to 30%, in addition to a possible 10% amplitude calibration uncertainty. A few case 
have >100% uncertainty in the column densities due to large uncertainties in the optical 
depth. 

We also computed spectra averaged within a 12" radius (roughly the beamwidth at of the 
IRAM 30m at 93 GHz) centered on the position of each protostar. We used these spectra to 
derive average abundances and column densities toward the protostars for comparison with 



previous single-pointing studies (e.g. Roberts & Millar 2007; Emprechtinger et al. 2009) 



3. Observational Results 

We obtained simultaneous N 2 D + and N 2 H + mapping of twelve protostellar envelopes, 
as well as interferometric data in both species for LI 157. The best example of spatially 
resolved N 2 D + and N 2 H + depletion is found toward L1157 in the interferometric data. We 
will discuss this source in detail first before presenting the rest of the sample. 



3.1. L1157 



The SMA continuum and N 2 D+ ( J = 3 -»■ 2) data (4.3"x7.1" resolution) for L1157 are 
shown in Figure 1, overlaid on the Spitzer IRAC 8 /im image. The 1.3 mm continuum data 
show that the source is strongly dominated by the point source emission coincident with the 
protostar, arising from the dense inner protostellar envelope. There is also extension along 
the outflow direction, likely due to warm gas associated with the inner outflow cavity. The 



outflow is nearly within the plane of the sky (Gueth et al. 1997) having an inclination of 
~85°. At lower intensity levels, the continuum emission is found to be extended along the 



dense envelope that was detected in extinction at 8 /xm (Looney et al. 2007 Tobin et al. 



2010) and in N 2 H + /NH 3 molecular line emission dChiang et al.||2010[ [Tobin et al.|201l] ). This 
extension along the extinction lane was previously seen in the continuum with single-dish 



bolometer mapping at 850 fim and 1.3 mm (Gueth et al. 2003), but not in interferometric 



continuum data (Gueth et al. 1997, 2003 Chiang et al. 2010 2012). The 3 mm continuum 



image from the PdBI is centrally peaked with no indication of extension along the envelope 



(Figure 1), consistent with Chiang et al. (2010). The detection of emission extended along 



the envelope with the SMA likely owes to the high-sensitivity of the observation, the short 
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uv-spacings in subcompact configuration, and the brighter dust continuum emission at 1.3 
mm as compared to 3 mm. 

The N2D + ( J = 3 — > 2) integrated intensity map (Figure 1), unlike the continuum data, 
is double-peaked, roughly symmetric about the protostar, with a local minimum coincident 
with the position of the protostar. The peaks are at radii of ~ 2000 AU (7") from the 
protostar and the emission closely correlates with opaque regions of the 8 /im extinction 
feature. In contrast, the N2H + integrated intensity map from the PdBI+30m has emission 
peaks that are located at R~1000 AU (3.3") from the protostar (Figure 1); the N2D + emission 
peaks are 1000 AU further out. The double-peaked N2H + emission appears at the same 



location in both the PdBI+30m image and the PdBI-only image from (Tobin et al. 2011). 
The symmetry about the protostar in N 2 D + and N 2 H + is striking, a strong indication that 
the protostar is affecting the abundance distribution of these molecules. 

We also conducted single-dish mapping of N2H + in the ( J = 3 — > 2) and ( J = 1 — > 0) 
transitions and N 2 D + in the (J = 3 — > 2) and (J = 2 — > 1) transitions. The integrated 
intensity maps of these lines are shown in Figure [2] and spectra of all four lines toward the 
location of the protostar are shown in Figure |3j The N2H + maps both appear centrally 
peaked and the ( J — 3 — > 2) transition appears slightly elongated along the envelope axis, 
consistent with the resolved PdBI data. Both N2D + maps show a double-peaked feature, 
coincident with the widely separated emission peaks observed in the SMA data. 

The single-dish data having centrally peaked N2H + (J = 3 — > 2) emission, double- 
peaked N2D + (J = 3 — > 2) and (J = 2 — > 1) rule-out any excitation anomaly, optical 
depth, or interferometric filtering being responsible for the apparent central depletion of 
N 2 D + emission around L1157. If the higher temperatures at smaller radii were causing 
the apparent central deficit of N 2 D+ (J = 3 -)> 2), then the N 2 H+ (J = 3 ->■ 2) should 
have appeared to be similarly double-peaked. Optical depth is ruled out because fitting the 
hyperfine transitions of N2D + finds them to be optically thin; the optically thin nature of 
the emission is apparent in the spectra shown in Figure 3 where the satellite lines emission 
lines are quite faint relative to the main lines. Therefore, we conclude that the deficit of 
N 2 D + emission toward the protostar and the offset between interferometric N 2 H + and N 2 D + 
peaks reflects a reduction of N2D + abundance relative to N2H + in the inner envelope. 

We used the interferometric N 2 H + and N 2 D + maps to construct column density maps 
of each species in order to construct a ratio map of N 2 D + to N 2 H + . The N 2 H + data have 
the single-dish data combined with the interferometer data to correct for missing flux. The 
N2D + data have resolved out some flux, apparent from the negative contours in Figure 1; 
however, the total interferometric flux is comparable to the single-dish flux indicating that 
most emission is recovered in the SMA map. The SMA data were sampled with the same 
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pixel scale and phase center as the PdBI data and we calculated the column densities for 
each transition following Goldsmith & Langer (1999). The excitation temperature for N 2 D + 
was taken from the hyperfine fits of N2H + (J = 1 — > 0) at each point because the N 2 D + 
( J = 3 — > 2) was optically thin and the excitation temperature could not be determined from 
hyperfine fitting. Figure [4] shows the SMA N2D + ( J = 3 — > 2) emission contours overlaid 
on the PdBI+30m N 2 H + (J = 1 — > 0) map and shows the column density ratio of N 2 D + 
relative to N 2 H + in the adjacent panel. We then take one dimensional radial cuts of the 
column densities and ratio and plot them in Figure |5j 

The line parameters for the average N 2 H + and N 2 D + spectra within a 12" radius around 
the protostar are given in Table 4 and the column densities and ratios are given in Table 5. 



We calculate that L1157 has 3 times the level of N 2 D + as measured by Emprechtinger et al. 



(2009), but a similar N 2 H + column density. This may result from us having mapped the 
area surrounding the protostar, while Emprechtinger et al. (2009) appear to have only taken 
spectra in single pointings centered on the protostar. In which case, the 16" beam at the 
frequency of N 2 D + ( J = 2 — > 1) would not have sampled the region of peak N 2 D + emission 
with only one pointing. The single-dish ratio is also an order of magnitude higher than the 
minimum observed in the interferometer data, mainly due to the N 2 H + column density being 
about and order of magnitude larger in the interferometer data. The lower column density 
in the single-dish data can be attributed to beam dilution and if we consider the column 
densities and ratios at 3600 AU radii (12"), then the interferometer data are consistent with 
the single-dish. 



3.2. Other Sources 

While we have given significant attention to L1157, observations of a larger sample 
show that N 2 D + and N 2 H + are often spatially differentiated and L1157 is not a special case. 
Figure 6 shows the N 2 H + emission maps of the full sample overlaid on the IRAC 8 /im 
images, compared with the N 2 D + emission maps. Figure 7 shows the column density maps 
of N 2 H + , N 2 D + , and the ratio of N 2 D + to N 2 H + for HH211 as an example and only the ratio 
maps are shown for the rest of the sample in Figure 8. The column densities are calculated 
by the method described in Section 2, using the integrated intensity corrected for optical 
depth and excitation temperature. The parameters of the N 2 H + and N 2 D + line fits are given 
in Table 4 and our column densities and ratios are given in Table 5. Our calculated column 



densities and ratios are consistent with those in Emprechtinger et al. (2009) within about 
a factor of two. The differences likely result from our exact centering on the protostar and 
our 12" radius for averaging the spectra to derive an overall column density and abundance 
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ratio. The sources are briefly discussed individually below. 



IRAS 03282+3035 - This source has a very dense envelope extended in the north- 
south direction, with nearly coincident N2H + and N2D + emission peaks. The ratio map 
indicates that the peak N2D + abundance is offset from the protostar and N2H + peak on one 
side only and to the south by 8". This offset seems real and is larger than the 2" receiver 
offset. 



HH211 - This protostar was the most heavily deuterated object in Emprechtinger et al. 



(2009 ) and this is true of our sample as well (see Table 5), with a comparable ratio of N2D + to 
N 2 H + . However, those measurements were toward the protostar position; the peak emission 
of N 2 H + and N 2 D + is offset to the southwest. Despite the N 2 H + emission peak being offset 
from the protostar, the column density peak is centered on the protostar (Figure 7). The 
N2D + column density peak is, however, offset southwest and the peak level of N2D + /N2H + is 
0.28 at (-9", 9"). 

L1521F - This protostar was classified as a Very Low Luminosity Object (VeLLO) by 



Bourke et al. (2006) and was previously regarded as a starless core (Crapsi et al. 2004). The 



N2D + emission peak may be slightly offset to the north relative to the N2H + ; the overall 
spatial distribution of N 2 D + appears quite flat; these maps are consistent with those of 
Crapsi et al. (2004). The column density ratio maps indicates that the N 2 D + abundance 



may be reduced toward the central protostar. 

IRAS 04325+2402 - This protostar is a multiple Class I source in Taurus, but with 
a large, dense core in close proximity, approximately located at (10, 25) in Figure 6 and 
offset from the protostar by ~63" (Scholz et al. 2010 Hogerheijde & Sandell 2000). The 



peak N 2 H + emission is found to be north of the protostar, coincident with the nearby dense 
core seen in 8 /im extinction, with emission extending southward toward the protostars. The 
N2D + emission is all located north of the protostar, coincident with the N2H + peaks. The 
detection of N 2 D + toward the protostar listed in Table 4 is marginal. 



L1527 - L1527 is also in Taurus and Emprechtinger et al. (2009) did not detect N 2 D 



emission toward this source. However, we find that the N2D + emission is concentrated north 
of the protostar by ~70" (10000 AU) and north of the N 2 H+ peak by ~35" (5000 AU), see 
Figure 6. Despite the large offset of N 2 D + from the protostar, the submillimeter emission 
and 8 /xm extinction are also more extended on the north-side of the protostar and coincide 



with the N 2 D+ emission (Tobin et al. 2010; Chandler & Richer 2000). Roberts & Millar 



(2007) also detected N 2 D + toward this source in a lower resolution study with comparable 
column density. Despite the Class status of L1527, it is very depleted in N 2 D + relative to 
others in the sample. This could result from the wide outflow cavities increasing the overall 
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envelope temperature and driving down the deuteration. 

RN043 - This source was the most distant in the sample, associated with Orion. 
Nonetheless, we did find that the N2D + emission was spatially offset from the N2H + emission, 
with the N2D + mostly on the west side of the protostar. The column density and ratio map 
is quite noisy due to the faintness of N 2 D + ; moreover, the abundance of N 2 D + relative to 
N 2 H + may be overestimated due to beam dilution of the N 2 H + map, given that RN043 has 
the lowest N 2 H + column density in the sample. 

HH108 IRS - In the case of this protostar, the N 2 D + and N 2 H + emission are both 
centered to the southeast of the protostar. Unfortunately, the signal-to-noise of the N 2 D + 
was too poor to make a useful column density ratio map. 

HH108 MMS - This is the nearby neighbor of HH108 IRS, but is apparently less 
evolved and lower luminosity. Its ratio of N 2 D + to N 2 H + is about 4 times greater than 
HH108 IRS, consistent with the apparent youth of this protostar relative to HH108 IRS. The 
N 2 D + emission from his protostar might be offset, but the noisy map makes this indefinite. 

L483 - L483 is the strongest N 2 H + source in the sample and most optically thick. 
The N 2 D + emission peak is offset by ~20" west of the protostar, with the highest N 2 D + to 
N 2 H + ratio at the location of peak N 2 D + emission. However, the column density ratio map 
(Figure 8) indicates that there could be higher column density N 2 D + coincident with the 
highest column density of N 2 H + . 

L673 - The subregion of L673 mapped coincides with L673-SMM2 as identified by 



Visser et al. (2002), harboring a small clustering of young stars. The sources marked in 



Figure 6 and 8 correspond to the Class 0/1 sources 13, 27, and 28 in Tsitali et al. (2010). 
The N 2 D + emission is clearly not peaked at the same location as the N 2 H + emission and 
the column density is substantially reduced. The average N 2 H + to N 2 D + ratios are fairly 
constant and low for the three identified sources. However, the northern most source (28) 
may be at the location of higher N 2 D + abundance as shown in the column density and ratio 
maps, but this is uncertain. 

L1152 - LI 152 is comprised of two cores seen in N 2 D + and N 2 H + emission linked by 
an apparent filament of material. The protostellar core in the south has strong N 2 D + and 
N 2 H + emission, but offset from the protostar. The peak N 2 D + to N 2 H + ratios in the south 
(protostellar) and north (starless) are 0.09 and 0.11 respectively. 

L1157 - We performed the same analysis on L1157 as the rest of the sample for com- 
pleteness, the drop in N 2 D + abundance is not apparent in the single-dish ratio map (Figure 
8), but there is a slight depression of emission in the N 2 D + (J = 2 — y 1) map. The aver- 
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age N 2 D + to N 2 H + ratio of 0.12 is consistent with averaging the minimum and maximum 
abundance ratio from the interferometer data. 

L1165 - Similar to L483, the N 2 D + emission is also shifted along the outflow in LI 165. 
The emission does seem to extend back toward the N 2 H + peak that is southeast of the 
protostar. The small area of N 2 D + emission makes a ratio map uninteresting, but the source 
has one of the lowest ratios of N 2 D + to N 2 H + from the single-point analysis; however this 
could be affected by beam dilution. 

Including L1157, eight of the thirteen sources show spatial shifts in the peak emission 
or column density of N 2 D + relative to N 2 H + . Although L673 is included, it could be argued 
that it is different since its environment is more complex due to having several young stars 
within a 0.1 pc diameter region. It is also important to note that we are limited by resolution 
in many cases, so some sources could have shifts that are not detectable in the 30m data. 
The N 2 H + and N 2 D + shifts can be regarded with high confidence because the data were 
taken simultaneously and the receivers are aligned to within 2" accuracy. The relative 
offset between N 2 H + and N 2 D + peaks are tabulated in Table 5. In addition to N 2 D + being 
shifted relative to N 2 H + , there are also offsets between the peak emission of N 2 H + and the 
protostars in six of the thirteen sources (HH108 IRS, L1152, L1165, HH211, L1527, and 
IRAS 04325+2402). 

While spatial differentiation between the species is quite prevalent, the location of the 
peaks relative to the protostar is often not as well-ordered as L1157. We observe a trend in 
the sample in which the apparently youngest sources, as classified from their cold bolometric 
temperatures {T M ) (Table 1) (HH211, HH108MMS, IRAS 03282+3035, and L1521F) do not 
show as strong spatial differentiation between N 2 D + and N 2 H + . The more evolved protostars 



(with warmer T bo i and wide outflows (Arce & Sargent 2006)) have the largest spatial offset 



between N 2 D+ and N 2 H+ (L1527 and IRAS 04325+2402). 

We also compare the N 2 D + to N 2 H + ratios to the bolometric luminosity and temperature 
in Figure [9] There is not a clear trend in bolometric luminosity, only that more luminous 
sources tend to have less deuteration while lower luminosity sources may be more highly 
deuterated or not. There does appear to be a trend of decreasing N 2 D + to N 2 H + ratio with 
increasing bolometric temperature, indicating that the deuteration may be decreasing with 
protostellar evolution. However, there is significant scatter at low bolometric temperatures. 
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3.3. Kinematics of N 2 D + versus N 2 H 



The difference in relative spatial distributions between N 2 D + and N 2 H + in most sys- 
tems is an indication that these molecules are not tracing the same material throughout 
the envelope. This suggests that the kinematic structure probed by these tracers is likely 



different. We compare the kinematics of the two tracers for L1157 in Figures 10 & [TTJ 
using both single-dish and interferometric kinematics respectively. For L1157, we notice a 
similar large-scale velocity gradient in both N 2 H + and N 2 D + ; however, near the protostar, 
where the N 2 D + abundance drops, there are clear differences. The line-of-sight velocities 
measured from N 2 D + on the west side of the envelope appear systematically larger by about 
0.1 km s _1 as compared to the N 2 H + data; this is present in both the single-dish and inter- 
ferometric data. This region is located near the blueshifted outflow cavity and could be due 



to N 2 H + being entrained by the outflow (Tobin et al. 2011), skewing the line profile toward 



blue-shifted velocities. Since the N 2 D + abundance is peaked at larger radii, it may not be 
affected by the outflow as strongly as N 2 H + . Thus, the N 2 D + could be tracing an increased 



velocity gradient due to rotation and/or infall (Tobin et al. 2011, 2012) that is masked by 
the outflow interaction with N 2 H + . 

The linewidth of N 2 H + is also found to be larger than N 2 D + within 10" of the protostar 
in both the single-dish and interferometer data. Thus, the signature of the outflow-envelope 
interaction is not as pronounced in N 2 D + due to its lower abundance on scales less than 2000 
AU. This is a clear indication that on small-scales, the N 2 D + is not tracing the same gas as 
N 2 H + . Moreover, this demonstrates that N 2 H + is able to probe closer to the protostar and 
disk forming region of the envelope during this stage of evolution. 



4. Chemical Model Comparison 

L1157 represents an ideal case to test chemical models for a protostellar source due to 
the higher resolution data, excellent signal-to- noise, and well-known orientation in the plane 



of the sky. Moreover, the relative symmetry of the filamentary envelope (Tobin et al. 2010), 
N 2 H + , and N 2 D + emission on scales less than 5000 AU make chemical modeling of this 
source more tractable. We will first review the critical chemical processes at work, affecting 
the formation and abundance of N 2 H + and N 2 D + . 
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4.1. Relevant Chemical Processes 



Deuterium chemistry in the cold (T~ 10 K) interstellar medium begins with the forma- 
tion of H 2 D + , dominated by the reaction 

H+ + HD ^ H 2 D + + H 2 + AE (1) 



where AE ~ 230 K (Herbst 1982). At low temperatures (T < 20 K) the forward reaction 



will dominate causing H 2 D + to have an abundance approaching that of Ht in cold, dense 



cores (Langer 1985). The high abundance of reactive H 2 D + enables subsequent deuterated 
molecules further down the reaction chain to have abundances elevated above the cosmic 

_ 5 i v w » 1 1 1 /" 1 1 v wv T~\ _ T— T - ! - nnr] 1 



^ Casein et al. 


2008 


Pagani et al. 


2009 



[D/H] ~ 10" 

dance of H 2 D + is limited by recombination with electrons and reactions with other molecules. 
The back reaction can also become non-negligible at T ~ 20 K, leading to decreased deuter- 



ation even before CO has been liberated from the dust (Emprechtinger et al. 2009). 



At T < 20 K, most CO is frozen-out onto dust grains at the densities of protostellar 
cores (n~10 5 cm~ 3 ), but once the protostar raises the temperature to > 20 K, the CO ice 
sublimates back to the gas phase. When CO is in the gas phase, all the available carbon not 
locked into dust grains will react, effectively shutting off the deuterium chemistry as initiated 
by H 2 D + ; an additional route to deuterium fractionation is available through CH 2 D + at T 
> 30 K, but this does not cause subsequent N 2 D + formation. Furthermore, the H 2 ortho- 
to-para ratio of the cloud can affect the overall deuterium fractionation; ortho-H 2 is in the 
J=l state with an energy above the ground of 170 K. This has enough energy to lower the 
potential barrier of the back reaction of Equation 1 even at low temperatures, making the 
overall level of deuteration sensitive to the H 2 ortho-to-para ratio ( |Flower et aL]|2006t |Pa gam 



et al. 2009) 



The primary formation route for N 2 H + is the gas phase reaction 



N 2 + H+ -> N 2 H H 



(2) 



and N 2 D + is understood to form via the same process, but with H 2 D + as the reactant. Thus, 
the ability to form N 2 D + is constrained by the conditions in which H 2 D + exists. Both N 2 D + 
and N 2 H + are destroyed by reactions with CO, therefore, they are only expected to be present 
with significant abundance in regions where CO has depleted. In low- mass/luminosity pro- 
tostellar systems, the temperatures only rise above 20 K at radii less than 1000 AU. Thus, 
the results from chemical modeling lead us to expect N 2 D + and N 2 H + to trace the same gas 
where CO is frozen-out onto dust grains (e.g. Emprechtinger et al.||2009 Pagani et al.| [2009). 



- 15 - 



4.2. Chemical Model 



We compare our data with the evolutionary chemical model of Lee et al. (2004). Since 



L1157 has formed a protostar, a model that accounts for the infalling material as well as the 
evolution of radiative heating is important. This model endeavors to simulate the chemistry 
of an isolated protostellar cloud through the prestellar phase as approximated by a Bonnor- 
Ebert sphere and then transition to the protostellar phase using the inside-out collapse model 
( Shu|1977 ) with an initially isothermal temperature of 10 K. The protostellar phase includes 



luminosity evolution following the prescription of (Young & Evans 2005) with radiative 



transfer calculated by DUSTY, but not including episodic accretion bursts. The chemistry 
of the infalling gas parcels is calculated discretely and the inward radial motion is tracked 
within the model. 

Since the chemical model calculates a ID radial abundance profile, we took several steps 
to make a more realistic observational comparison. We used the input radial density profile 
from the chemical model and the calculated abundance profile to construct a circular enve- 
lope midplane. We then summed the density along the line of sight through the envelope 
to calculate the predicted column density distribution in order to compare with the obser- 



vations. Since L1157 is likely a filament rather than a sheet (Tobin et al. 2010), we also 



tapered the circular sheet by a Gaussian along the assumed line of sight with scale height 
parameters of 500, 750, 1000, and 2000 AU in order to simulate a filamentary density distri- 
bution. However, we note that the chemical abundances calculated for a ID density profile 
will not exactly map to a filamentary envelope, but should be a reasonable approximation. 
Our goal was not to fit L1157 in detail with the model, but to look for a general agreement 
with the observed radial distribution and abundance ratio of the two species. 



In Figure 12 we show a plot of the N2D + to N2H + column density ratios as a function 
of radius from the model. We clearly see the drop in abundance of the two species at radii 
less than 1000 AU due to reactions with CO and their abundances peak at intermediate 
radii between 1000 and 5000 AU, with decreasing abundances at larger radii. The model 
that best fits the abundance ratio of L1157 has evolved for 300,000 yr after collapse and the 
luminosity of the model (5 L Q ) is comparable to what is observed. While the abundance 
ratio of N 2 D + relative to N 2 H + matches the data reasonably well, the column density and 
abundance peaks of N 2 H + and N 2 D + are nearly coincident with each other, inconsistent with 



the data shown in Figure 12 Earlier and later timesteps of the model also show coincident 
peaks between these two molecules. Thus, predicted ratio of N 2 D + and N 2 H + is consistent 
with the model, but the absolute abundances/column densities are inconsistent. 

It is important to point out that the column density peak of N 2 D + is in agreement 
with the observed N 2 D + column density peak and the N 2 H + column density peak of the 
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model is what disagrees with the observations. The inability of the model to reproduce the 
abundance peak offsets between the N 2 D + and N 2 H + may indicate that the models may not 
be capturing some physical process; we will discuss this further in the following section. 



Discussion 



It has been previously shown that the depletion of N2D + relative to N2H + is expected 
in protostars (Roberts & Millar 2007 Emprechtinger et al. 2009), but the ~1000 AU spa- 



tial differentiation between the peaks of the two species is not matched by our chemical 



modeling, nor that of Emprechtinger et al. (2009). Moreover, the only other study to recog- 



nize the spatial offset of N2H + and N2D + was Friesen et al. (2010) with observations of the 



protocluster Oph B2. However, such spatial offsets have not been observed toward individ- 
ual/isolated protostellar systems previously. Spatially resolved observations of deuterated 
molecules relative to their normal counterparts offer crucial tests to chemical models and 
further our understanding of the evolution of deuterated molecules in protostellar systems. 



5.1. N 2 H+ Depletion Zones 



As pointed out in the introduction, a N 2 H + depletion zone had been observed previously 



toward the Class protostar IRAM 04191 (Belloche & Andre 2004). This source is recognized 



as a Very Low Luminosity Object or VeLLO, having an internal luminosity of < 0.1 L . 
The depletion zone in IRAM 04191 has a radius of 1600 AU, larger than that of L1157 
and other protostars observed by Chen et al. (2007) and Tobin et al. (2011). However, 



the disappearance of N 2 H + could be explained by either freeze-out onto dust grains, as 



apparently observed in the pre-stellar core B68 (Bergin et al. 2002), or the destruction of 



N 2 H + due to CO being released from dust grains as a result of protostellar heating (e.g. Lee 
et al.|200~4 ). Such a large depletion zone in IRAM 04191 produced via N2H + destruction from 



evaporated CO would not work given the low system luminosity and Belloche & Andre (2004) 



suggested freeze-out of N 2 onto dust grains as the likely mechanism for N 2 H + depletion. 



Alternatively, Lee (2007) constructed an episodic accretion model that could also explain 



the size of the depletion zone with a recent burst of luminosity in IRAM 04191. L1157 and 
the other protostars known to have resolved N2H + depletion zones are of order 5 to 10 L Q . 
These protostars could evaporate CO from dust grains out to sufficiently large radii with 
equilibrium heating to explain the size of the depletion zones. 
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5.2. Cause of N 2 H+ and N 2 D+ Emission Offsets 



The inability of the chemical model used here (Lee et al. 2004) and others in the litera- 



ture ( Emprechtinger et al.|2009 ) to reproduce the relative offset of the column density peaks 
for both species is troubling at first. However, this may not mean that our understanding of 
the chemical processes forming and destroying the molecules is wrong. The presence of the 
central protostar, raising the gas and dust temperatures above ~10 K introduces significant 
complexities to the chemistry within ~2000 AU of the protostar, as compared to the pre- 
stellar phase. Moreover, since we are considering the radial abundances of two molecules, 
we must determine whether it is the N 2 D + or N 2 H + that is problematic. 

As outlined in the previous section, the reaction that forms H 2 D + , leading to the for- 
mation of N 2 D + , is only active at temperatures less than 20 K; at greater temperatures 
the H 2 D + reacts with H 2 and is converted back to Hg" and HD (Equation 1). Furthermore, 
T ~ 20 K is also the temperature at which CO is evaporated off of dust grains and can 
freely react with and destroy N 2 D + and N 2 H + . A C 18 map of L1157 detects an emission 



peak at the protostar and very little CO at radii larger than 1200 AU (J0rgensen et al. 



2007), consistent with the location of the N 2 H + abundance peaks. CO will destroy N 2 D + 
as quickly as it destroys N 2 H + , so reactions with CO will affect the abundance distribution 
of both molecules in the same manner. In order to have an abundance peak of N 2 H + inside 
that of N 2 D + , one needs to have a range of radii where the H 2 D + abundance is reduced and 



CO is not significantly present in the gas phase. Figure [13] shows a schematic diagram of 
the spatial zones in L1157 where N 2 H + , N 2 D + , and the two molecules governing the spatial 
distribution of emission are located (CO and H 2 D + ). Of these molecules, only observations 
of H 2 D + are lacking. 

The relevant processes in the protostellar envelope that could affect the spatial distribu- 



tion of molecules are protostellar heating, outflows, and infall. Emprechtinger et al. (2009) 



presented chemical models that considered protostellar heating and the H 2 D + abundances 
becoming reduced at radii somewhat larger than CO sublimation; however, the shifts of 
N 2 D + peak abundances were only ~200 AU outside of the N 2 H + peak and unable to explain 
the 1000 AU offset in LI 157. Replenishment of N 2 H + and N 2 D + via infall from the outer 



envelope is built into Lee et al. (2004) model, but N 2 H + and N 2 D + remain coincident at all 



time steps. Tobin et al. (2011) showed that the protostellar outflow is entraining N 2 H + at 
the 1000 AU scale; however, such entrainment would affect N 2 D + in the same way that it 
affects N 2 H + . Moreover, N 2 D + does not appear to be significantly affected by the outflow 



as evidenced by its narrower linewidth than N 2 H + (Figure 10). 



Another possibility considered was episodic accretion (e.g. Kenyon et al.|1990 Dunham 
et al. 2010). If there was an accretion outburst, it could temporarily heat the envelope to 
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T > 20 K out to ~2000 AU. Then all the H 2 D+ could have been converted back to H^ 
and HD. Moreover, CO would be liberated from the dust grains destroying the N 2 H + and 



N2D + present (Lee 2007 Visser & Bergin 2012). After the luminosity decreases and CO 



freezes back out, N2H + would be able to form directly, while N2D + would have to wait until 
H2D + has reformed. Thus, we considered the possibility that N2H + could reform faster 
than N 2 D + since it is a one-step process and N 2 D + takes two steps to form. While this 



scenario is attractive, work by Lee (2007); Visser & Bergin (2012) shows that the chemistry 
is governed by the CO freeze-out timescale after an accretion outburst and N2D + and N2H 4 " 
should reform in lockstep. Thus, episodic accretion alone cannot adequately account for the 
observed spatial distribution of N2H + and N2D + . 

An additional scenario we considered is an H 2 ortho-to-para (OPR) gradient in the 
envelope. The steady-state OPR in cold (T ~ 10 K) clouds is calculated to be ~2.7xl0~ 3 
Flower et al. (2006), whereas H2 forms with OPR = 3 from statistical equilibrium. OPR 



values larger than steady-state would decrease the deuteration by enabling the back reaction 
of H 2 D + formation to be active (Equation 1). Since the protostar likely went through a cold 
starless phase with a duration of 10 4 to 10 6 yr the H 2 OPR throughout the cloud was likely 



around the equilibrium value before collapse (Parise et al. 2011). Once the protostar forms 



and begins to heat the inner envelope, the OPR ratio could increase; Sipila et al. (2010) 



shows a sharp increase at T > 15 K. Such an OPR ratio gradient going into the warm inner 
envelope could explain the large region of N 2 D + depletion. However, it is unclear if this 
could reproduce the observations given that the OPR is still quite low (~10 -4 ). There is 
also a question of timescales given that the steady-state timescale in Sipila et al. ( |2010 ) is 
~10 5 yr, longer than the infall timescale for material to fall in from 2000 AU to 1000 AU (~ 
10 4 yr assuming a 0.5 M protostar). While there are some problems with this scenario, it 
cannot be ruled out. 

A final scenario to produce the observed abundance offsets is a CO evaporation temper- 
ature is higher than 20 K. The temperatures for H2D + destruction and CO evaporation are 
quite close and if the CO evaporation temperature was even just a few degrees higher than 



20 K, this could be enough to explain the observed emission offsets, as was true in the Em- 



prechtinger et al. (2009) models. We illustrate this in Figure 13, where if the CO remained 



mostly frozen-out between 20 to 25 K, this could explain the emission offsets by enabling 
N2H + to be present in this warmer region at smaller radii than N2D + . It is possible that the 
CO evaporation temperature could vary since the physical state of CO ice can be much more 
complex. The CO may be mixed with other ices: NH 3 , C02,CH 3 OH, H 2 0, etc. The binding 



energy of CO to other ices tends to be strong than CO to CO (Collings et al. 2004 Fayolle 



et al. 2011), therefore ice mixtures having slightly larger evaporation temperatures could 



explain the observations. However, it is uncertain how mixed the ice is because some species 
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tend to freeze-out more quickly, forming ice layers rather than a single mixed layer (Oberg 



et al. 2011). Aside from the caveat of unknown ice mixtures, a slightly altered evaporation 



temperature of CO ice due to mixing would explain our observations quite simply without 
needing episodic accretion and/or H 2 OPR gradients. But, episodic accretion could cause 
more ice mixing and conversion of CO ice to CO2 may build an ice mixture with a large 
enough binding energy to remain frozen at T > 20 K (Kim et al. 2012). 

In summary, the abundance peak offsets between N2D + and N2H + could be explained by 
either an H2 OPR gradient in the inner envelope or an increased CO evaporation temperature 
due to ice mixtures. The chemical model of Lee et al. (2004) does not calculate a temperature 
dependent H 2 OPR, which may explain why it does not reproduce the abundance peak offset 
between N 2 D + and N 2 H + . At the same time, the model also does not allow for ice mixtures 
to increase the evaporation temperature CO, which could cause the N 2 H + peak to move 
inward. However, the chemical model does have its N 2 D + column density /abundance peak 
at roughly the same radius as the observations, this offers some evidence for the increased 
CO evaporation temperature scenario, i.e. the deuterium chemistry produced the expected 
results for N 2 D + but N 2 H + is divergent. Though, we cannot conclusively distinguish between 
these two scenarios as both seem likely to occur and could be simultaneously contributing 
to make the shift in abundance peaks stronger. 



5.3. Envelope and N 2 D + /N 2 D + Emission Morphology 

An additional effect to consider is the viewing angle and three dimensional envelope 
density structure in relation to the observed N 2 D + emission morphology. In the case of LI 157, 



there is evidence that the envelope is filamentary rather than sheet-like (Tobin et al. 2010). 
A filamentary envelope can have two effects on the apparent N 2 H + and N 2 D + depletion 
zones. First, the depletion zones are likely more apparent in a filamentary geometry than in 
a axisymmetric sheet, due to less emission superimposed along the line of sight in than an 



axisymmetric object; this is demonstrated in the column density plots in Figure 12 Also, 
the orientation of a filament within the plane of the sky can change the apparent size of 
the depletion regions; however, viewing angle effects would only make the depletion regions 
smaller and not larger. Thus, no complicated viewing angle effects could give the illusion of 
large N 2 D + depletion zones outside those of N 2 H + . In the case of LI 157, we cannot attribute 
the observed morphology to observational bias due to the symmetry of the depletion zones 
and abundance peaks. 



The rest of the sample is not as well-ordered as LI 157 in terms of envelope struc- 
ture and/or emission morphology. L1527, L1165, L1152, HH108 IRS, HH211, and IRAS 
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04325+2402 have N 2 H + emission peaks that are not coincident with the protostar, but the 
N 2 D + emission is never peaked closer to the protostar that N 2 H + . L483 and L1165 show 
N 2 D + emission peaked in the direction of the outflow rather than perpendicular to it. In 
these cases of complex emission morphology, the heating of the protostar and the envelope 
structure are both contributing to the complex emission morphology. L483 shows N 2 H + , 



NH 3 , and 850 /im emission that are oriented at a 45° angle to the outflow direction (Fuller 



& Wootten 2000 Shirley et al. 2000), roughly the same direction as the N 2 D + emission. 
Therefore, the N 2 D + emission apparently extended along the outflow is simply shifted in the 
direction of the high-density envelope structure and the highest density material is simply 
not perpendicular to the outflow. 

Moreover, the envelopes with less complex structure near the protostar in 8 /mi extinc- 
tion and N 2 H+ (i.e. L1521F, IRAS 03282+3035, HH211, and L1157) do not have irregular 
N 2 D + distributions. These sources, except for L1157, appear to have N 2 D + emission peaks 
possibly shifted to one side of the envelope. Thus, we can conclude that if the envelope has 
highly irregular structure, then we expect this to create an irregular distribution of N 2 D + . 
Furthermore, the frequent offsets of N 2 D + emission from the protostars, even though it may 
be coincident with the N 2 H + peak, means that N 2 D + is not tracing the densest material in 
protostellar cores. Thus, caution must be exercised when interpreting the kinematics of the 
deuterated species without observations of both as deuterated species cannot be relied upon 
to trace the densest and/or same material as non-deuterated in the presence of temperature 
gradients. 



5.4. Chemical Evolution 



It is important to stress that even though N 2 D + and N 2 H + are not always spatially 
coincident as previously assumed and indicated by chemical models, our results do not 
disagree with the generalized picture of chemical evolution in protostellar cores (e.g. Lee 



et al. 2004 Roberts k Millar 2007 Emprechtinger et al. 2009). Specifically, the reduced 



abundance of N 2 D + relative to N 2 H + as protostar evolve is consistent with our data as 
shown in Figure [9] and the previous work; moreover, the spatial offsets of N 2 D + relative 
to N 2 H + appear to be a product of this chemical evolution and contribute to the relative 
decrease of N 2 D + relative to N 2 H + . Given the two scenarios we outline for the cause of 
the peak offsets, further chemical modeling that accounts for infall, H 2 OPR gradients, and 
changes in the CO evaporation temperature is needed to fully understand the exact cause. 
Sensitive, spatially-resolved observations of both C 18 and H 2 D + could determine if the 
increased evaporation temperature is likely, as we would expect to see a gap between the 
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two species if this were the case. 



6. Summary and Conclusions 

We have presented a detailed study of the relative spatial distribution of N2D + and 
N2H" 1 " in Class 0/1 protostellar envelopes from interferometric observations of LI 157 and 
single-dish data of thirteen sources (L1157 included). We find that N2H + and N2D + often 
do not have the same column density distribution, given that the emission peaks of the two 
species are not spatially coincident, contrary to predictions by previous chemical models. 

1. The depletion regions of both N 2 H+ (J = 1 -> 0) and N 2 D + (J = 3 ->■ 2) in the 
envelope around the Class protostar LI 157 are well-resolved in interferometric observations 
with the PdBI and SMA respectively. Both N 2 D + and N 2 H + are double-peaked with the 
N 2 D+ peaks at R~2000 AU, outside the peaks of N 2 H+ which are at R~1000 AU. Follow-up 
single-dish observations from the IRAM 30m of N 2 H+ (J = 3 -> 2), N 2 D+ (J = 3 -»• 2), and 
( J — 2 — >• 1) confirm the spatial offsets of the species. 

2. Subsequent observations of additional protostellar envelopes in N 2 H + ( J = 1 — > 0) 
and N 2 D + (J = 2 — > 1) with the IRAM 30m show frequent offsets between the emission 
peaks of N 2 D + and N 2 H + in eight of thirteen sources (L1157 included). Some sources without 
offsets could be due to a lack of spatial resolution; moreover, the additional sources with 
spatial offsets are not as well-ordered as LI 157. 

3. The large (~1000 AU) offset in abundance peaks of N 2 H + and N 2 D + is not matched 



by the dynamic chemical model, nor previously published static models from (Emprechtinger 



et al. 2009). In such models, the abundance ratio of N2D + and N2H + falls due to heating 



by the protostar, but the abundance peaks remain closely coincident. 

4. We suggest that the emission offsets are best explained by either an increased CO 
evaporation temperature due to ice mixtures or an H2 ortho-to-para ratio gradient. In the 
case of an increased CO evaporation temperature, there would be a region where H2D + is 
not present in the gas phase, but CO is still frozen-out, enabling N 2 H + have its emission 
peak be inside of the N 2 D + peak. However, an H 2 ortho-to-para ratio gradient would appear 
observationally similar and both effects may be happening. 
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Fig. 1. — L1157- Data from the SMA and PdBI overlaid on the IRAC 8 //m image (inverse 
grayscale). The top left panel shows the SMA 1.3 mm continuum emission (contours) which 
peaks on the protostar but is also extended at low surface brightness levels along the envelope 
viewed in 8 /zm extinction; contours are ±2cr, ±3, 6, 9, 12, 22, 30, 50, 150cr, where a 
= 2.16 mJy beam -1 . The top right panel shows the SMA N2D + (J = 3 — > 2) integrated 
intensity emission, the contours start at and increase in units of ±3<r, where a = 0.07 K 
km s -1 . The N 2 D + emission shows a local minimum on the protostar position and peaks at 
R ~ 2000 AU. The bottom left panel shows the PdBI 3.4 mm continuum emission (contours) 
which peaks on the protostar; contours are ±2cr, ±3, 6, 12, 24, 36, 96a, where a = 0.24 
mJy beam -1 . The bottom right panel show the PdBI + IRAM 30m N2H + (J = 1 — > 0) 
integrated intensity emission; the contours start at and increase in units of ±30<r, where 
a = 0.042 K km s -1 . The N 2 H + emission also shows a local minimum at the location of 
the protostar, but the peak emission is located at R~1000 AU and inside the peak N2D 4 " 
emission. The SMA beam is 7.1" x 4.3" and the PdBI beam is 3.4" x 3.3". 



-27- 



N a H + (1-0) Integrated Intensity K-km/s N 2 H + (3-3) Integrated Intensity K-km/s 




Aa ["] Act ["] 



Fig. 2. — LI 157- Single-dish mapping from the IRAM 30m. The top panels show the 
N2H + ( J — 1 — y 0) and (J = 3 -> 2) emission, respectively. The bottom panels show the 
N 2 D + (J = 2 — >■ 1) and (J = 3 — >• 2) emission. The N2H + emission is centrally peaked 
and the N 2 H + (J = 3 — > 2) appears elongated normal to the outflow, consistent with the 
double-peaked N 2 H + (J = 1 — > 0) map from the PdBI (Figure 1). The N 2 D + emission then 
appears double-peaked in both transitions, consistent with the SMA N 2 D + observations. 
The contours and beam in each panel are drawn as follows: N 2 H + ( J = 1 — > 0) lOer, 30cr, 
50(7, and rr=0.027 K km s" 1 , 27" beam; N 2 H+ (J = 3 -> 2) 4rj, 6a, and (7=0.14 K 
km s" 1 , 9" beam; N 2 D+ (J = 2 ->■ 1) 2rr, 4a, and a=0.05 K km s" 1 , 16" beam; N 2 D+ 
(J = 3 -» 2) 2a, 4a, and a=0.12 K km s" 1 , 11" beam. 
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Fig. 3. — Example spectra of N 2 H + and N 2 D + taken toward L1157. The spectra are averaged 
over an aperture radius of 12" centered on the protostar position. Faint hyperfine lines are 
visible in the N 2 H + (J = 3 — >■ 2) at and 5 km s^ 1 and in N 2 D+ (J = 2 ->■ 1) at -3 and 
6 km s -1 . The signal-to- noise is too low to see hyperfine lines in N 2 D + (J = 3 — >■ 2). The 
negative features in the N 2 H + (J = 3 — >■ 2) at -5 and -12 km s -1 are frequency switching 
artifacts. 
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Fig. 4. — Comparison of N 2 H + and N 2 D + emission and column density ratio map. The left 
panel shows the N2D + integrated intensity (contours) overlaid on the PdBI N2H + integrated 
intensity (color scale). The right panel shows the column density ratio of N2D 4 " to N2H + . 
The color stretch in the left panel goes from to 21 K km s _1 and the contours start 3a and 
increase with this interval (a = 0.07 K km s _1 ). 
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Fig. 5. — N 2 H + and N 2 D + column density and column density ratio plots. The left panel 
shows the column density of N 2 H + and N 2 D + derived from the integrated intensity maps. 
The right panel shows the N 2 D + /N 2 H + column density ratio. Note the clear 1000 AU 
separation of N 2 H + and N 2 D + column density peaks. The error bars on the N 2 H + data are 
the combination of noise, error in excitation temperature, and error in optical depth. The 
error bars on the N 2 D + data only include noise since the emission is optically thin and we 
adopt the excitation temperature of the N 2 H + data. The solid lines are taken from the west 
side of the envelope and the dotted lines are from the east side. 
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Fig. 6. — Single-dish integrated intensity maps of N 2 H + ( J = 1 — » 0) and N2D + ( J = 2 — > 1) 
for the protostars identified in each panel. The left panels show the N 2 H + (J = 1 — > 0) 
emission overlaid on the IRAC 8 /im image, the middle panels show the N 2 H + ( J = 1 — > 0) 
as contours and grayscale, and the right panels show the N 2 D + (J = 2 — > 1) emission as 
grayscale and contours. Eight of thirteen well-resolved sources show spatial differentiation 
between the N 2 H + and N 2 D + emission peaks; however, they are not as regular as in LI 157. 
All N 2 H + contours start at 10a and increase by 10cr and all the N 2 D + contours start at 3a 
and increase in 3a intervals; the values for a are given in Table 3 for each source. The IRAM 
30m beam for the N 2 H+ ( J = 1 0) data is 27" and 16" for the N 2 D+ ( J = 2 ->• 1) data. 
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Fig. 7. — Spatial distribution of N 2 H + column density (left panel), N 2 D + (center panel), 
and their ratios (right panel) for HH211. Notice the column density peak of N 2 H + at the 
location of the protostar, while the emission peak was offset in Figure 6. 
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Fig. 8. — Column density ratio maps for all sources except HH211 (Figure 7), L1165, and 
HH108. 
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Fig. 9. — N 2 D + to N 2 H + ratio plotted versus bolometric luminosity (left panel) and bolomet- 
ric temperature (right panel). Higher luminosity sources tend to have a lower N 2 D + /N 2 H + 
ratio, while lower luminosity sources have a range of deuteration. There appears to be a 
trend in bolometric temperature versus N 2 D + /N 2 H + ratio with deuteration decreasing with 
evolution, consistent with (Emprechtinger et al. 2009). 
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Fig. 10. — Kinematic comparison of N2H + and N2D + from the IRAM 30m data of LI 157. 
The top row shows the IRAM 30m N 2 H + (J = 1 — > 0) integrated intensity (left panel), 
line-center velocity (middle panel), and linewidth (right panel). The bottom row shows the 
same information as the top but for N 2 D + (J = 2 — > 1). The line-center velocities for N 2 D + 
are about 0.1 km s -1 more redshifted on the west side of the envelope than on the east side 
and the linewidths of N 2 D + are lower than the N 2 H + linewidths toward the protostar. 



-38- 




Fig. 11. — Kinematic comparison of interferometric N2H + and N2D + data of L1157. The top 
row shows the PdBI N2H + (J = 1 — > 0) integrated intensity (left panel), line-center velocity 
(middle panel), and linewidth (right panel). The bottom row shows the same information as 
the top but for SMA N 2 D + (J = 3 — > 2). The velocity structure of N 2 H + in the line-center 
velocities and linewidth is much more complex than N 2 D + due to the outflow interaction 
with the envelope (Tobin et al. 2011). The more redshifted gas on the west side of the 
envelope also appears in the N 2 D + (J = 3 — > 2) which is also present in the single-dish 
N 2 D+. 
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Fig. 12. — Abundance ratio and integrated column densities of N 2 H + and N 2 D + compared 
to models. The abundance ratio of the ID model and tapered sheet (filament) resemble the 
abundance ratio of the data; however, the absolute column densities do not match the data. 
The left panel shows the abundance ratio of N 2 D + to N 2 H + from the ID chemical model 
calculation (dashed-line), integrated through a circular sheet with a radius of 15000 AU 
(triple dot-dashed line) and an approximated filamentary envelope, produced by tapering 
the circular sheet by a Gaussian with a = 500 AU (single dot-dashed line). The data are 
shown as the solid and dotted lines from the east and west sides of the envelope respectively. 
The right panel shows the column density of N 2 H + and N 2 D + derived from integrated maps, 
assuming that N 2 D + has the same excitation temperature as N 2 H + . The data are shown as 
the solid and dotted lines from the east and west sides of the envelope respectively and the 
model data are shown as the dashed or dot-dashed lines. Note the coincidence of column 
density peaks in the models for N 2 D + ; moreover, the peaks are most apparent in the filament 
model rather than the sheet. 
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Fig. 13. — Schematic diagram for spatial distribution of N 2 H + and N 2 D + in L1157. The 
upper panel shows the suggested relative abundances of each species as a function of radius 
and temperature and the bottom panel identifies the temperature range and which species 
are present in the gas phase for the case of increased CO evaporation temperature. Inside 
of 1000 AU, CO is in the gas phase, evaporated off the dust grains, and destroys the N 2 H + . 
At radii between 1000 AU and 2000 AU the CO may remain frozen on the dust grains at 
the expected temperatures of 20 to 25 K; this enables N 2 H + to be abundant, but H 2 D + is 
destroyed due to the back reaction operating (Equation 1) at these temperatures, keeping 
N 2 D + abundances low. At radii greater than 2000 AU, the temperature is below 20 K and 
H 2 D + will be present, enabling N 2 D + to form and have its abundance peaks at this radius. 
If there is an H 2 OPR gradient rather than an increased CO evaporation temperature, N 2 D + 
and H 2 D + may only be present at T < 15 K and CO could then evaporate at T > 20 K and 
still show the same abundance offset. 
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Table 2. Observed Species 



Transition 


Frequency 
(GHz) 


Observatory 


Beam width 
(") 


F e// 


B eff 


N 2 H+ (J = 1 ->• 0) 


93.1737637 


IRAM 30m 


26 


0.95 


0.81 






PdBI D, C-array 


3.4 x 3.3 






N 2 H+ (J = 3 -> 2) 


279.511832 


IRAM 30m 


9 


0.88 


0.53 


N 2 D+ (J = 2 ->• 1) 


154.217206 


IRAM 30m 


16 


0.93 


0.74 


N 2 D+ (J = 3 ->■ 2) 


231.321966 


IRAM 30m 


10 


0.94 


0.63 






SMA sub-compact 


7.1 x 4.3 







Note. - The single-dish data were rescaled to main beam brightness tempera- 
tures (T m (,) where T m (, = F e yyT^/B e yy. F e // is the forward efficiency, account- 
ing for losses by the secondary support structure, and B e ff is the main beam ef- 
ficiency, the ratio of the main beam solid angle over the entire antenna pattern 
(http://www.iram.es/IRAMES/mainWiki/Iram30mEfficiencies). 



Table 3. IRAM 30m Observations 



Source 


RA 


Dec 


Date 


T t 

J- sys 


Channel 


Widtht 






(J2000) 


(J2000) 


(UT) 


(K) 


(kHz) 


(K km s- 1 ) 


IRAS 03282+3035 


03:31:19.8 


+30:45:37 


25 Oct 2009 


115, 130 


20, 


40 


0.04, 0.035 


HH211-mm 


03:43:55.3 


+32:00:46 


24 Oct 2009 


110, 120 


20, 


40 


0.05, 0.04 


L1521F 


04:28:39.0 


+26:51:37 


24 Oct 2009 


102, 110 


20, 


40 


0.04, 0.025 


IRAS04325+2402 


04:35:36.6 


+24:08:54 


25 Oct 2009 


105, 120 


20, 


40 


0.04, 0.035 


L1527 


04:39:54.0 


+26:03:22 


26 Oct 2009 


100, 130 


20, 


40 


0.03, 0.015 


RN043 


05:32:18.2 


+12:50:01 


24 Oct 2009 


115, 135 


20, 


40 


0.04, 0.03 


L483 


18:17:30.3 


-04:39:40.6 


25 Oct 2009 


135, 180 


20, 


40 


0.06, 0.05 


HH108 


18:35:42.7 


-00:33:07 


23 Oct 2009 


140, 230 


20, 


40 


0.05, 0.06 


L673-SMM2 


19:20:25.4 


+11:19:44 


23/24 Oct 2009 


140, 210 


20, 


40 


0.04, 0.04 


L1152 


20:35:47.17 


+67:53:22 


24 Oct 2009 


120, 150 


20, 


40 


0.04, 0.035 


L1157 


20:39:05.6 


+68:02:13.2 


25 Oct 2009 


110, 115, 320, 290 


20, 20, 


40, 40 


0.027, 0.05, 0.12, 0.14 


L1165 


22:06:50.8 


+59:03:06.5 


25 Oct 2009 


110, 120 


20, 


40 


0.04, 0.035 



Note. — Observations of protostcllar envelopes taken with the IRAM 30m telescope. The positions reflect the map center and not 
necessarily the protostar position. 

^System temperatures, channel widths, and intensity uncertainties are given for 93 GHz and 154 GHz respectively. For L1157, values 
for 231 GHz and 279 GHz are also given. 
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Table 5. Column densities of N?H + and N?D + Ratios and Offsets 



Source 



N(N 2 H+) 
(cm" 2 ) 



N(N 2 D+) 
(cm" 2 ) 



N(N 2 D+)/N(N 2 H+) N2H+/N2D+ offset 

(") 



IRAS03282 


2 


66 


X 


10 13 


± 


3 


25 


X 


10 12 


3.42 


X 


10 12 


± 


6 


.52 


X 


10 11 


0.129 


± 





.029 


8 


HH211 


3. 


.41 


X 


10 13 


± 


4. 


17 


X 


10 12 


7.57 


X 


10 12 


± 


1 


.58 


X 


10 12 


0.222 


± 





.054 


20 1 " 


L1521F 


7. 


42 


X 


10 13 


± 


1 


.72 


X 


10 12 


2.09 


X 


10 12 


± 


2. 


00 


X 


10 11 


0.028 


± 





.003 





IRAS04325 


1. 


51 


X 


10 13 


± 


3 


.97 


X 


10 12 


< 1.' 


5 x 10 10 










< 0.001 








IRAS04325N 


2. 


83 


X 


10 13 


± 


5 


42 


X 


10 12 


4.05 


X 


10 11 


± 


4, 


.89 


X 


10 10 


0.014 


± 





.003 





L1527 


1. 


59 


X 


10 13 


± 


3 


.53 


X 


10 12 


< 1.' 


5 x 10 10 










< 0.001 






35 


L1527N 


1. 


56 


X 


10 13 


± 


3 


.78 


X 


10 12 


1.85 


X 


10 11 


± 


3. 


59 


X 


10 10 


0.012 


± 





.004 




RN043 


3 


18 


X 


10 12 


± 


1 


35 


X 


10 12 


1.79 


X 


10 11 


± 


7, 


02 


X 


10 10 


0.056 


± 





.033 


16 


L483 


1 


40 


X 


10 14 


± 


6 


.17 


X 


10 12 


6.92 


X 


10 12 


± 


1. 


60 


X 


10 12 


0.049 


± 





.012 


16 


HH108 IRS 


2. 


01 


X 


10 13 


± 


5 


.41 


X 


10 12 


4.03 


X 


10 11 


± 


8. 


27 


X 


10 10 


0.020 


± 





.007 


Q a 


HH108 MMS 


5, 


99 


X 


10 12 


± 


2 


55 


X 


10 12 


5.32 


X 


10 11 


± 


7, 


23 


X 


10 11 


0.089 


± 





.127 




L673S 


1 


89 


X 


10 13 


± 


5. 


.05 


X 


10 12 


3.47 


X 


10 11 


± 


5. 


89 


X 


10 10 


0.018 


± 





.006 


18 


L673M 


1 


77 


X 


10 13 


± 


4 


.56 


X 


10 12 


2.54 


X 


10 11 


± 


4. 


84 


X 


10 10 


0.014 


± 





.005 




L673N 


1 


.57 


X 


10 13 


± 


4. 


05 


X 


10 12 


2.60 


X 


10 11 


± 


3 


.98 


X 


10 11 


0.017 


± 





.026 




L1152 


2. 


09 


X 


10 13 


± 


5. 


20 


X 


10 12 


1.49 


X 


10 12 


± 


8, 


79 


X 


10 11 


0.071 


± 





.046 





L1157 


1 


82 


X 


10 13 


± 


4 


.49 


X 


10 12 


2.19 


X 


10 12 


± 


5. 


.58 


X 


10 11 


0.120 


± 





.043 


3.5 


L1165 


3 


76 


X 


10 12 


± 


1 


.97 


X 


10 12 


3.75 


X 


10 10 


± 


2. 


03 


X 


10 10 


0.010 


± 





.008 


8 



tTaken from the column density maps, emission maps are spatially coincident. 

a Shift between protostar and N 2 H+ is observed even though N 2 H+ and N 2 D+ are spatially coincident. 



